##########################################################################################

library(dplyr)
library(ggplot2)
library(data.table)
library(RColorBrewer)
library(optparse)

##########################################################################################

option_list <- list(
    make_option(c("--gene"), type = "character") ,
    make_option(c("--mut_rate_gene_file"), type = "character") ,
    make_option(c("--mut_rate_point_file"), type = "character") ,
    make_option(c("--images_path"), type = "character")
)

if(1!=1){
    
    gene <- "MUC6"
    work_dir <- "~/20220915_gastric_multiple/dna_combine/"
    mut_rate_gene_file <- paste(work_dir,"/images/mutRate/MutRate.tsv",sep="")
    mut_rate_point_file <- paste(work_dir,"/images/mutRate/MutRate.RecurrentPoint.tsv",sep="")
	images_path <- paste(work_dir,"/images/mutRatePlot",sep="")

}

###########################################################################################

parseobj <- OptionParser(option_list=option_list, usage = "usage: Rscript %prog [options]")
opt <- parse_args(parseobj)
print(opt)

gene <- opt$gene
mut_rate_gene_file <- opt$mut_rate_gene_file
mut_rate_point_file <- opt$mut_rate_point_file
images_path <- opt$images_path

dir.create(images_path , recursive = T)

###########################################################################################

col <- c(
  brewer.pal(9,"YlGnBu")[6],
  rgb(234,106,79,alpha=255,maxColorValue=255),
  rgb(203,24,30,alpha=255,maxColorValue=255),
  rgb(255,0,0,alpha=255,maxColorValue=255)
  )

names(col) <- c("IM" , "IGC" , "DGC" , "GC")

###########################################################################################

dat_mutRateGene <- data.frame(fread( mut_rate_gene_file ))
dat_mutRatePoint <- data.frame(fread( mut_rate_point_file ))

###########################################################################################

dat_mutRateGene <- subset(dat_mutRateGene , Hugo_Symbol==gene)
dat_mutRatePoint <- subset(dat_mutRatePoint , Hugo_Symbol==gene)
dat_mutRatePoint$Hugo_Symbol <- dat_mutRatePoint$vid
dat_mutRatePoint <- dat_mutRatePoint[,-3]

dat_plot <- rbind( dat_mutRateGene , dat_mutRatePoint )
dat_plot$Class <- factor( dat_plot$Class , levels = c("IM" , "GC" , "IGC" , "DGC") , order = T )
dat_plot$value_text <- paste0( round(dat_plot$MutRate , 2) * 100 , "%") 


###########################################################################################
trans <- function(num){
    up <- floor(log10(num))
    down <- round(num*10^(-up),2)
    text <- paste("p == ",down," %*% 10","^",up)
    return(text)
}

###########################################################################################
## 计算P值
dat_plot$p.value = ""
dat_plot$p_text = ""

for(geneN in unique(dat_plot$Hugo_Symbol)){

    print(geneN)

    tmp_1 <- subset( dat_plot , Hugo_Symbol == geneN & Class %in% c("GC") )
    tmp_2 <- subset( dat_plot , Hugo_Symbol == geneN & Class %in% c("IM") )

    if(nrow(tmp_1)==0){
        tmp_1 <- tmp_2
        tmp_1$MutNum <- 0
        tmp_1$MutRate <- 0
        tmp_1$value_text <- ""
    }

    if(nrow(tmp_2)==0){
        tmp_2 <- tmp_1
        tmp_2$MutNum <- 0
        tmp_2$MutRate <- 0
        tmp_2$value_text <- ""
    }

    tmp <- rbind( tmp_1 , tmp_2 )

    tmp_fisher <- matrix(c(tmp$MutNum , tmp$SampleNum - tmp$MutNum) , ncol = 2)
    p <- fisher.test(tmp_fisher)$p.value

    if( p < 0.001 ){
        p_text <- trans(p)
    }else{
        p_text <- paste0( "p == " , round(as.numeric(p) , 3) ) 
    }

    dat_plot[dat_plot$Hugo_Symbol == geneN & dat_plot$Class %in% c("GC" , "IM"),"p.value"] <- p
    dat_plot[dat_plot$Hugo_Symbol == geneN & dat_plot$Class %in% c("GC" , "IM"),"p_text"] <- p_text
}

dat_plot$Hugo_Symbol <- factor( dat_plot$Hugo_Symbol , 
    levels = unique(dat_plot[order(dat_plot$MutRate , decreasing=T),"Hugo_Symbol"]) , order = T)


###########################################################################################

plot <- ggplot( data = dat_plot , aes( x = Class , y = MutRate , fill = Class ))+
geom_bar(position = "stack", stat = "identity") + 
theme_bw()+
labs(x="",y="Mutation Rate")+
facet_grid(.~Hugo_Symbol) +
theme(panel.grid = element_blank())+
scale_fill_manual(values=col) +
ylim(0,1)+
geom_text(aes(label=p_text , y = 0.9 ,x = 2.5),parse = TRUE,size=4)+
geom_text(aes(label=value_text) , position=position_stack(vjust = 0.5) , size=2.5 , color="black")+
theme(panel.background = element_blank(),#设置背影为白色#清除网格线
            legend.position ='right',
            legend.title = element_blank() ,
            panel.grid.major=element_line(colour=NA),
            plot.title = element_text(size = 12,color="black",face='bold'),
            legend.text = element_text(size = 12,color="black",face='bold'),
            axis.text.y = element_text(size = 12,color="black",face='bold'),
            axis.title.x = element_text(size = 12,color="black",face='bold'),
            axis.title.y = element_text(size = 12,color="black",face='bold'),
            strip.text.x = element_text(size = 7 , face = 'bold'),
            axis.ticks.x = element_blank(),
            axis.text.x = element_text(size = 8,color="black",face='bold') ,
            axis.line = element_line(size = 0.5))

out_name <- paste0( images_path , "/MutRate_" , gene , ".pdf" )
ggsave( out_name , plot , width = (3 + length(unique(dat_plot$Hugo_Symbol))) , height = 5 )